R Notebook

Import Libraries

install.packages("astsa") 

Error in install.packages : Updating loaded packages

library(astsa)

Examine Cmort dataset

Average weekly cardiovascular mortality in Los Angeles County; 508 six-day smoothed averages obtained by filtering daily values over the 10 year period 1970-1979.

• plot cmort dataset to take a look at the data
• plot lag1 of cmort dataset

Simperasmata apo ta Dedomena 1. Paratiroume oti to dataset exei mia elafros fthinousa tasi 2. Paratiroume oti iparxei mia epoxikotita sta dedomena (lamvanontas ypoopsin oti kathe xronos exei 52 vdomades elenxontas tis 501 vdomades pou exoun kataxorithei sto dataset)

Ektimisi pithanon montelwn Apota PACF kai ACF plots paratiroume ta eksis: 1. To ACF grafima fenete na fthinei (oxi omos ekthetika) 2. Sto PACF grafima vlepoume 2 spikes gia lag 1 kai 2

# load c mort dataset and make sure that it is a time series 
x=ts(cmort)

# plot cmort data
plot.ts(x, main="Cmort Dataset")

# plot cmort lag 1
lag1.plot(cmort,1)

# plot ACF for  24 lags
acf(x, xlim=c(1,24))

# plot PACF or  24 lags
acf2(x)

     [,1] [,2] [,3] [,4]  [,5]  [,6]  [,7]  [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  0.77 0.77 0.68 0.65  0.58  0.53  0.48  0.41 0.39  0.32  0.28  0.23  0.18  0.13  0.11  0.04  0.01 -0.03 -0.07 -0.08 -0.10 -0.12 -0.13 -0.13
PACF 0.77 0.44 0.03 0.03 -0.01 -0.05 -0.02 -0.05 0.05 -0.08 -0.03  0.00 -0.06 -0.06  0.05 -0.10 -0.04 -0.01 -0.03  0.03  0.04 -0.02  0.00  0.01
     [,25] [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33]
ACF  -0.15 -0.16 -0.14 -0.17 -0.14 -0.15 -0.11 -0.10 -0.07
PACF  0.00 -0.04  0.07 -0.08  0.03  0.01  0.07  0.03  0.03

Differecing Dataset

• Dedomenou oti paratirithike mia mikri fthinousa tasi tha xriastei na paroume tis protes diafores gia tin eksalipsi autis
• Episis exei paratirithei kai mia epoxikotita ana peripou 52 periodous opote oi protes diafores mporei na voithisoun kai stin eksalipsi autis.

Symperasmata Meta tin efarmogi ton proton diaforwn apo ta diagramata ACF kai PACF prokiptei oti ena AR(1) montelo tha itan katalilo gia tin epilisi tou provlimatos

xdiff=diff(x,1)
# plot cmort data
plot.ts(xdiff, main="Cmort Dataset After First Diffs")

plot(xdiff, type="b")

# plot ACF for  24 lags
acf(xdiff)

# plot PACF or  24 lags
acf2(xdiff)

      [,1]  [,2]  [,3]  [,4]  [,5]  [,6] [,7]  [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  -0.51  0.21 -0.13  0.06 -0.02 -0.01 0.02 -0.09 0.12 -0.10  0.04  0.00  0.01 -0.07  0.12 -0.09  0.02 -0.01 -0.06  0.01  0.01 -0.03 -0.02  0.03
PACF -0.51 -0.07 -0.07 -0.03  0.01 -0.01 0.02 -0.09 0.04 -0.01 -0.04  0.02  0.02 -0.09  0.06  0.00 -0.04 -0.01 -0.08 -0.08 -0.01 -0.04 -0.05 -0.03
     [,25] [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33]
ACF  -0.01 -0.07  0.13 -0.14  0.08 -0.08  0.04 -0.04  0.03
PACF  0.01 -0.10  0.05 -0.07 -0.04 -0.10 -0.06 -0.05 -0.03

Erotima 1)

Έστω xt η χρονοσειρά (cmort). Προσαρμόστε τα δεδομένα σε ΑR(2) μοντέλο χρησιμοποιώντας γραμμική παλινδρόμηση. Στη συνέχεια με το προσαρμοσμένο μοντέλο που βρήκατε να κάνετε πρόβλεψη σε ορίζοντα 4 βδομάδων και να βρείτε ένα 95% διάστημα εμπιστοσύνης για την πρόβλεψη.

xlag1 = lag(xdiff,-1)
xlag2 = lag(xdiff,-2)
# create new dataset
y = cbind(xdiff,xlag1,xlag2)
# drop NA
y<-y[complete.cases(y),]
k<-data.frame(y)

xdiff <- k$xdiff
xlag1 <- k$xlag1
xlag2 <- k$xlag2

fit = lm(xdiff~xlag1+xlag2, data = data.frame(y))
summary(fit)

Call:
lm(formula = xdiff ~ xlag1 + xlag2, data = data.frame(y))

Residuals:
     Min       1Q   Median       3Q      Max 
-19.1228  -3.8222  -0.1499   3.6291  20.5285 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.03631    0.25885  -0.140    0.888    
xlag1       -0.53798    0.04448 -12.096   <2e-16 ***
xlag2       -0.06652    0.04455  -1.493    0.136    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.817 on 502 degrees of freedom
Multiple R-squared:  0.2587,    Adjusted R-squared:  0.2557 
F-statistic: 87.58 on 2 and 502 DF,  p-value: < 2.2e-16

confint(fit, level=0.95)

                 2.5 %      97.5 %
(Intercept) -0.5448849  0.47225714
xlag1       -0.6253609 -0.45059497
xlag2       -0.1540441  0.02100059

Pragmatopoiooume provlepsi gia 4 vdomades (diladi 4 fores) xrisimopoiontas ta apotelesmata tis proigoumenis provlepsis san input stin epomeni

results <- c()
for (i in 1:4) {
  
  # 334  -1.89 -12.77   8.36
  # 335 x0 -1.89 -12.77
  # 335 x1 -12.77 x0
  # 335 x2 x0 x1
  # 335 x3 x1 x2
 last_row = tail(y,1)
 newData =  data.frame(xlag1 = last_row[1], xlag2 = last_row[2])

 pred_val = predict(fit, newdata = newData)
 results <- append(results, pred_val)

 # update y array
 y <- rbind(y, c(pred_val, last_row[1], last_row[2]))
}

print("Forcast results for 4 weaks")

[1] "Forcast results for 4 weaks"

print(results)

          1           1           1           1 
 1.35756663 -0.50455898  0.14481998 -0.08065968 

Symperasmata

Paratiroume oti pernoume ta idia apotelesmata me tin methodo sarima.for pou parexei to package astasa - sarima.for results: 1.36352984 -0.50021224 0.14907469 -0.07605003 - our model results: 1.35756663 -0.50455898 0.14481998 -0.08065968

# dimiourgia AR(1) model
sarima (xdiff, 2, 0, 0)

initial  value 1.908985 
iter   2 value 1.801760
iter   3 value 1.763215
iter   4 value 1.759680
iter   5 value 1.759665
iter   6 value 1.759665
iter   6 value 1.759665
iter   6 value 1.759665
final  value 1.759665 
converged
initial  value 1.758285 
iter   2 value 1.758284
iter   2 value 1.758284
iter   2 value 1.758284
final  value 1.758284 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
          ar1      ar2    xmean
      -0.5379  -0.0666  -0.0182
s.e.   0.0444   0.0444   0.1610

sigma^2 estimated as 33.65:  log likelihood = -1604.5,  aic = 3216.99

$degrees_of_freedom
[1] 502

$ttable
      Estimate     SE  t.value p.value
ar1    -0.5379 0.0444 -12.1198  0.0000
ar2    -0.0666 0.0444  -1.4981  0.1347
xmean  -0.0182 0.1610  -0.1129  0.9101

$AIC
[1] 6.370287

$AICc
[1] 6.370382

$BIC
[1] 6.403749

sarima.for(xdiff, 4, 2, 0, 0)

$pred
Time Series:
Start = 506 
End = 509 
Frequency = 1 
[1]  1.36352984 -0.50021224  0.14907469 -0.07605003

$se
Time Series:
Start = 506 
End = 509 
Frequency = 1 
[1] 5.800736 6.586576 6.712074 6.729726

Erotima 2)

Δημιουργείστε n=500 παρατηρήσεις από ένα ARMA μοντέλο με τον τύπο Xt= 0.9Xt-1 + Wt - 0.9Wt-1 οπου Wt~iidN(0,1)

Paratirisis 1. Arxika paratiroume oti ta ACF kai PACF dixnoun oti ta dedomena den mporoun na anaparastathoun mono me ena AR i ena MA montelo kati pou mas odigei sto simperasma oti ena ARMA montelo einai i lisi 2. episis fenete oti ta residuals tou ACF den einai statistika simantika (kati pou ipodilonei oti to montelo mas einai kalo) 3. ta p-values fenete na einai ola mi statistika simantika (kati pou ipodilonei oti to montelo mas einai kalo) 4. kai to Q-Q plot dixnei oti i provlepsis tou montelou einai kales.

Ta parapano simperasmata einai anamenomena dioti ta dedomena paraxthikan simfona me to sigkekrimeno montelo.

# Genetate 500 observations using arima.sim
gen_data <- arima.sim(model = list(ar = 0.9, ma = 0.9), n = 500)
# make sure that gen_data is time series 
gen_data=ts(gen_data)

# plot series
plot.ts(gen_data, main="Generated Series")

# plot ACF and PACF 
acf2(gen_data)

     [,1]  [,2] [,3]  [,4] [,5]  [,6] [,7]  [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  0.96  0.87 0.79  0.73 0.68  0.64 0.60  0.56 0.53  0.52  0.50  0.48  0.46  0.43  0.41   0.4  0.39  0.40  0.40  0.40  0.40  0.41  0.40  0.39
PACF 0.96 -0.49 0.35 -0.14 0.16 -0.11 0.07 -0.02 0.17 -0.08  0.08 -0.11  0.00  0.03  0.09   0.0  0.11 -0.02  0.07 -0.04  0.07 -0.02 -0.03  0.00
     [,25] [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33]
ACF   0.38  0.37  0.36  0.35  0.33  0.32  0.30  0.28  0.26
PACF  0.01 -0.02  0.06 -0.09  0.01  0.02 -0.04 -0.07  0.07

# apply data on ARIMA(1,0,1)
sarima (gen_data, 1, 0, 1)

initial  value 1.495208 
iter   2 value 0.036900
iter   3 value 0.012950
iter   4 value 0.012270
iter   5 value 0.012092
iter   6 value 0.011868
iter   7 value 0.011860
iter   8 value 0.011859
iter   9 value 0.011859
iter  10 value 0.011857
iter  11 value 0.011847
iter  12 value 0.011844
iter  13 value 0.011844
iter  14 value 0.011844
iter  14 value 0.011844
iter  14 value 0.011844
final  value 0.011844 
converged
initial  value 0.012885 
iter   2 value 0.012835
iter   3 value 0.012805
iter   4 value 0.012803
iter   5 value 0.012794
iter   6 value 0.012793
iter   7 value 0.012776
iter   8 value 0.012765
iter   9 value 0.012764
iter  10 value 0.012762
iter  11 value 0.012760
iter  12 value 0.012760
iter  13 value 0.012760
iter  14 value 0.012760
iter  15 value 0.012760
iter  16 value 0.012760
iter  17 value 0.012759
iter  18 value 0.012759
iter  19 value 0.012759
iter  20 value 0.012759
iter  21 value 0.012759
iter  22 value 0.012759
iter  23 value 0.012759
iter  24 value 0.012759
iter  25 value 0.012759
iter  26 value 0.012759
iter  27 value 0.012758
iter  28 value 0.012758
iter  29 value 0.012758
iter  30 value 0.012758
iter  31 value 0.012758
iter  32 value 0.012758
iter  33 value 0.012758
iter  34 value 0.012758
iter  34 value 0.012758
iter  34 value 0.012758
final  value 0.012758 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
         ar1     ma1    xmean
      0.9071  0.8873  -0.2743
s.e.  0.0187  0.0303   0.8979

sigma^2 estimated as 1.017:  log likelihood = -715.85,  aic = 1439.7

$degrees_of_freedom
[1] 497

$ttable
      Estimate     SE t.value p.value
ar1     0.9071 0.0187 48.5421  0.0000
ma1     0.8873 0.0303 29.3226  0.0000
xmean  -0.2743 0.8979 -0.3055  0.7601

$AIC
[1] 2.879393

$AICc
[1] 2.87949

$BIC
[1] 2.91311

Erotima 3)

Θεωρείστε τα δεδομένα oil Προσαρμόστε τα δεδομένα με ένα ARIMA (p,d,q) κατάλληλο μοντέλο και κάντε το διαγνωστικό έλεγχο του μοντέλου καθώς και την αξιολόγηση αυτού

Dataset Info Annual oil production (millions of tonnes), Saudi Arabia, 1965-2010.

Symperasmata 1. Apo to diagrama tis xronosiras paratiroume oti exoume anodiki tasi. Kati pou mas odigei sto na efarmosoume protes diafores gia na apalipsoume tin tasi stin xronosira. 2. episis kai to ACF diagramma to epiveveonei auto kathos den fthinei ekthetika.

# load oil dataset and make sure that it is a time series 
x=ts(oil)

# plot oil data
plot.ts(x, main="Oil Dataset")

# plot oil lag 1
lag1.plot(oil,1)

# plot ACF and PACF 
acf2(x)

     [,1]  [,2]  [,3]  [,4] [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  0.99  0.99  0.98  0.97 0.95  0.94  0.93  0.92  0.90  0.89  0.87  0.86  0.84  0.82  0.81  0.79  0.77  0.76  0.74  0.72  0.71  0.69  0.68  0.66
PACF 0.99 -0.17 -0.04 -0.10 0.00 -0.04 -0.03 -0.01 -0.14 -0.03 -0.04 -0.02 -0.01  0.05 -0.05  0.00  0.06 -0.04 -0.04  0.02  0.09  0.05 -0.04  0.04
     [,25] [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34]
ACF   0.65  0.64  0.63  0.61  0.61  0.59  0.59  0.58  0.57  0.56
PACF  0.00  0.10 -0.02  0.11 -0.02 -0.06  0.00  0.03  0.04  0.00

1st Differencing

Tha dokimasoume 2 diaforetikes paralages gia differecing: Arxika dokimazoume tis protes diafores. Xrisimopoiontas tis opoies prokoiptoun ta eksis. - Apo ta ACF kai PACF den antiprosopeuoun kapoio AR i kapoio MA montelo. Opote kataligoume sto simerasma na dokimasoume to ARMA(1,1) pano sta dedomena pou prokiptoun apo tis 1es diafores. - Apo tin aksiologisi tou parapano montelou ta apotelesmata den einai poli entharintika dioti paratiroume statistika simantika reciduals tou lag 7, 25 kai 25 - Episis kai to Ljung-Box den einai katholou entharintiko afou ola ta p-values einai mi statistika simantika.

xdiff=diff(x,1)
# plot oil data
plot.ts(xdiff, main="Oil Dataset After First Diffs")

# plot ACF and PACF for first diff data
acf2(xdiff)

     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24] [,25]
ACF  0.18 0.06 0.11 0.05 0.08 0.04 0.03 0.19 0.09  0.07  0.08  0.04 -0.04  0.05  0.00 -0.07  0.02  0.06 -0.04 -0.08 -0.08  0.03 -0.07 -0.08 -0.13
PACF 0.18 0.03 0.10 0.02 0.06 0.01 0.01 0.18 0.02  0.04  0.02  0.01 -0.09  0.05 -0.03 -0.10  0.03  0.04 -0.07 -0.09 -0.03  0.05 -0.07 -0.01 -0.12
     [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34]
ACF  -0.02 -0.15 -0.09  0.01 -0.05 -0.06 -0.07 -0.06 -0.01
PACF  0.02 -0.13  0.01  0.05 -0.02  0.00 -0.05  0.02  0.02

# apply data on ARIMA(1,0,1)
sarima (xdiff, 1, 0, 1)

initial  value 0.952571 
iter   2 value 0.942681
iter   3 value 0.937165
iter   4 value 0.937152
iter   5 value 0.937039
iter   6 value 0.936783
iter   7 value 0.934747
iter   8 value 0.933927
iter   9 value 0.933726
iter  10 value 0.933677
iter  11 value 0.933316
iter  12 value 0.932849
iter  13 value 0.932728
iter  14 value 0.931401
iter  15 value 0.931353
iter  16 value 0.931266
iter  17 value 0.931094
iter  18 value 0.930978
iter  19 value 0.930956
iter  20 value 0.930950
iter  21 value 0.930944
iter  22 value 0.930941
iter  23 value 0.930934
iter  24 value 0.930933
iter  25 value 0.930932
iter  26 value 0.930930
iter  27 value 0.930929
iter  28 value 0.930929
iter  28 value 0.930929
final  value 0.930929 
converged
initial  value 0.930069 
iter   2 value 0.930068
iter   3 value 0.930068
iter   4 value 0.930068
iter   5 value 0.930068
iter   6 value 0.930068
iter   6 value 0.930068
iter   6 value 0.930068
final  value 0.930068 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
         ar1      ma1   xmean
      0.8739  -0.7702  0.0629
s.e.  0.0600   0.0779  0.1968

sigma^2 estimated as 6.423:  log likelihood = -1277.86,  aic = 2563.72

$degrees_of_freedom
[1] 541

$ttable
      Estimate     SE t.value p.value
ar1     0.8739 0.0600 14.5683  0.0000
ma1    -0.7702 0.0779 -9.8897  0.0000
xmean   0.0629 0.1968  0.3198  0.7493

$AIC
[1] 4.712719

$AICc
[1] 4.712801

$BIC
[1] 4.744329

Log Differencing

Epidi ta dedomena mas fenete na akolouthoun ekthetiki tasi tha efarmosoume diafora ton logarithmon sto dataset oil. Meta tin diaforisi paratiroume: - Oti kai to ACF kai to PACF midenizonte meta to 4o lag kati mou mas ipodilonei tin iparksi enos ARMA(1,1) montelou gia ta dedomena pou prokipsan apo tin diaforisi. - efarmozontas to parapano montelo pernoume kalitera ACF residuals apo tis protes diafores kathos episis kai kalitero Ljung-Box - Episis kai to Q-Q plot dixnei na pragmatopoiei kales provlepsis akriopoiontas auto to montelo

Symperasmata: Apo ta parapano kataligoume oti gia to dataset Oil tha efarmosoume ena ARIMA montelo me: - AR taksis 1. - MA taksis 1. - Kai logarithmikes diafores.

# Calculate approximate oil returns
oil_returns <- diff(log(oil))

# Plot oil_returns. Notice the outliers.
plot(oil_returns)

# plot ACF and PACF for first diff data
acf2(oil_returns)

     [,1]  [,2] [,3]  [,4] [,5]  [,6]  [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  0.13 -0.07 0.13 -0.01 0.02 -0.03 -0.03 0.13 0.08  0.02  0.01     0 -0.02  0.06 -0.05 -0.09  0.03  0.05 -0.05 -0.07  0.04  0.09 -0.05 -0.08
PACF 0.13 -0.09 0.16 -0.06 0.05 -0.08  0.00 0.12 0.05  0.03 -0.02     0 -0.03  0.09 -0.07 -0.06  0.01  0.04 -0.05 -0.05  0.05  0.06 -0.06 -0.05
     [,25] [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34] [,35] [,36] [,37] [,38] [,39] [,40] [,41] [,42] [,43] [,44] [,45] [,46] [,47] [,48]
ACF  -0.07  0.00 -0.11 -0.07  0.02 -0.02 -0.03 -0.05 -0.03  0.00 -0.09 -0.01 -0.04 -0.01  0.02 -0.01 -0.06  0.01  0.00 -0.01  0.04  0.01  0.05  0.07
PACF -0.08  0.02 -0.11  0.01  0.00 -0.01 -0.05 -0.04  0.02  0.02 -0.08  0.02 -0.04  0.04 -0.01 -0.01 -0.05  0.03 -0.03  0.00  0.08  0.00  0.05  0.01
     [,49] [,50] [,51] [,52] [,53] [,54] [,55] [,56] [,57] [,58] [,59] [,60] [,61] [,62] [,63] [,64] [,65] [,66] [,67] [,68] [,69] [,70] [,71] [,72]
ACF  -0.01 -0.03  0.01 -0.04 -0.04 -0.03     0 -0.01 -0.10 -0.01 -0.05 -0.04 -0.03  0.01  0.01 -0.01 -0.04  0.02     0 -0.01 -0.03 -0.02 -0.05 -0.01
PACF  0.04 -0.08  0.01 -0.07  0.00 -0.06     0 -0.06 -0.11  0.01 -0.09 -0.01 -0.04  0.04 -0.01  0.00 -0.04  0.03     0  0.00 -0.04 -0.02 -0.04  0.00
     [,73] [,74] [,75] [,76] [,77] [,78] [,79] [,80] [,81] [,82] [,83] [,84] [,85] [,86] [,87] [,88] [,89] [,90] [,91] [,92] [,93] [,94] [,95] [,96]
ACF  -0.01 -0.02  0.01  0.02  0.04 -0.01  0.03  0.02 -0.04 -0.01  0.02  0.03  0.01  0.03  0.08 -0.04 -0.02  0.01 -0.04  0.05  0.07 -0.04  0.02  0.05
PACF -0.01  0.00  0.02 -0.01  0.04 -0.02  0.08 -0.03 -0.03 -0.03  0.03 -0.03 -0.02  0.03  0.04 -0.09 -0.01 -0.02 -0.03  0.03  0.05 -0.11  0.02 -0.01
     [,97] [,98] [,99] [,100] [,101] [,102] [,103] [,104] [,105] [,106] [,107] [,108] [,109] [,110] [,111] [,112] [,113] [,114] [,115] [,116] [,117]
ACF   0.01  0.00  0.01   0.04   0.01  -0.03  -0.04  -0.01   0.02   0.01   0.01   0.06   0.08   0.04   0.02   0.01   0.03   0.02  -0.02  -0.04  -0.01
PACF  0.02 -0.03  0.06   0.01  -0.05   0.02  -0.03   0.01   0.00   0.04  -0.01   0.07   0.04   0.04   0.00   0.05  -0.01   0.00  -0.04  -0.03  -0.03
     [,118] [,119] [,120] [,121] [,122] [,123] [,124] [,125] [,126] [,127] [,128] [,129] [,130] [,131] [,132] [,133] [,134] [,135] [,136] [,137]
ACF    0.04   0.05  -0.02  -0.02   0.03   0.01  -0.04  -0.08   0.02   0.00  -0.04   0.01   0.02   0.01   0.02   0.00  -0.01   0.00  -0.03  -0.06
PACF   0.02   0.04  -0.01  -0.04  -0.01   0.03  -0.03  -0.07   0.00  -0.02  -0.04   0.01   0.01  -0.01   0.02   0.05   0.02   0.01   0.02  -0.02
     [,138] [,139] [,140] [,141] [,142] [,143] [,144] [,145] [,146] [,147] [,148] [,149] [,150] [,151] [,152] [,153] [,154] [,155] [,156] [,157]
ACF    0.01  -0.02  -0.02   0.02  -0.01  -0.03   0.00   0.00  -0.04  -0.01  -0.02  -0.04  -0.04   0.01   0.01   0.04   0.03   0.01   0.05   0.01
PACF   0.04   0.01  -0.03  -0.02   0.02  -0.01   0.02  -0.01   0.02  -0.02  -0.03  -0.02  -0.01   0.02  -0.01   0.04   0.03  -0.04   0.03   0.00
     [,158] [,159] [,160] [,161] [,162] [,163] [,164] [,165] [,166] [,167] [,168] [,169] [,170] [,171] [,172] [,173] [,174] [,175] [,176] [,177]
ACF   -0.06   0.02   0.05  -0.02   0.05   0.00  -0.01      0  -0.01  -0.02  -0.01   0.00  -0.03  -0.01  -0.02  -0.02   0.04  -0.01  -0.03   0.02
PACF  -0.05   0.02   0.03   0.00   0.03   0.01   0.00      0   0.02   0.01  -0.01   0.03  -0.03   0.00  -0.04   0.00   0.02  -0.03  -0.01   0.01
     [,178] [,179] [,180] [,181] [,182] [,183] [,184] [,185] [,186] [,187] [,188] [,189] [,190] [,191] [,192] [,193] [,194] [,195] [,196] [,197]
ACF    0.01  -0.01  -0.01  -0.04   0.07  -0.01  -0.04   0.05  -0.02  -0.01   0.01  -0.05  -0.04  -0.01   0.01   0.04  -0.01   0.00   0.06  -0.06
PACF   0.02  -0.01   0.00  -0.04   0.08  -0.05   0.02  -0.01  -0.02   0.03   0.00  -0.03  -0.04   0.01  -0.01   0.07  -0.01   0.02  -0.01  -0.03
     [,198] [,199] [,200] [,201] [,202] [,203] [,204] [,205] [,206] [,207] [,208]
ACF   -0.02   0.02      0   0.00   0.00  -0.04   0.00   0.04   0.04   0.04   0.01
PACF   0.01   0.00      0  -0.02  -0.01  -0.05  -0.01   0.02   0.04   0.05  -0.01

# Assuming both P/ACF are tailing, fit a model to oil_returns
sarima(oil_returns, 1, 0, 1)

initial  value -3.057594 
iter   2 value -3.061420
iter   3 value -3.067360
iter   4 value -3.067479
iter   5 value -3.071834
iter   6 value -3.074359
iter   7 value -3.074843
iter   8 value -3.076656
iter   9 value -3.080467
iter  10 value -3.081546
iter  11 value -3.081603
iter  12 value -3.081615
iter  13 value -3.081642
iter  14 value -3.081643
iter  14 value -3.081643
iter  14 value -3.081643
final  value -3.081643 
converged
initial  value -3.082345 
iter   2 value -3.082345
iter   3 value -3.082346
iter   4 value -3.082346
iter   5 value -3.082346
iter   5 value -3.082346
iter   5 value -3.082346
final  value -3.082346 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
          ar1     ma1   xmean
      -0.5264  0.7146  0.0018
s.e.   0.0871  0.0683  0.0022

sigma^2 estimated as 0.002102:  log likelihood = 904.89,  aic = -1801.79

$degrees_of_freedom
[1] 541

$ttable
      Estimate     SE t.value p.value
ar1    -0.5264 0.0871 -6.0422  0.0000
ma1     0.7146 0.0683 10.4699  0.0000
xmean   0.0018 0.0022  0.7981  0.4252

$AIC
[1] -3.312109

$AICc
[1] -3.312027

$BIC
[1] -3.280499

Erotima 4)

Θεωρείστε τα δεδομένα globtemp. Προσαρμόστε τα δεδομένα με ένα ARIMA (p,d,q) κατάλληλο μοντέλο και κάντε το διαγνωστικό έλεγχο του μοντέλου καθώς και την αξιολόγηση αυτού. Στη συνέχεια κάντε πρόβλεψη για τα επόμενα 10 χρόνια

Analisi Dedomenon Global-mean monthly […] and annual means, 1880-present, updated through most recent month.

# load globtemp dataset and make sure that it is a time series 
x=ts(globtemp)

# plot oil data
plot.ts(globtemp, main="Global Temp Dataset")

# plot oil lag 1

Restarting R session...

Differecing

Apo tin parapano analisi ton dedomenon fenete oti ta dedomena exoun auksitiki tasi. Gia ton logo auto tha efarmosoume protes diafores prokimenoun na apalipsoume tin tasi auti. Efarmozontas ta ACF kai PACF sta dedomena afou paroume tis protes diafores gia na eksalipsoume tin tasi prokoiptoun ta eksis 3 pithana mpontela pou mporoume na xrisimopoieisoume: 1. ARIMA(1,1,1) auto to montelo dokimastike dedomenou oti kai to ACF kai to PACF pavoun kai ta 2 na einai statistika simantika kai den paratirite kapoio apo ta 2 na fthinei ekthetika 2. ARIMA(0,1,2) auto to montelo dokimastike dedomenou ACF grafima exei 2 spikes 3. ARIMA(3,1,0) auto to montelo dokimastike dedomenou PACF grafima exei 3 spikes

xdiff=diff(x,1)
# plot globtemp data
plot.ts(xdiff, main="Global Temp Dataset After First Diffs")

# plot ACF and PACF for first diff data
acf2(xdiff)

      [,1]  [,2]  [,3] [,4]  [,5]  [,6]  [,7] [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22]
ACF  -0.24 -0.19 -0.08 0.20 -0.15 -0.03  0.03 0.14 -0.16  0.11 -0.05  0.00 -0.13  0.14 -0.01 -0.08     0  0.19 -0.07  0.02 -0.02  0.08
PACF -0.24 -0.26 -0.23 0.06 -0.16 -0.09 -0.05 0.07 -0.09  0.11 -0.03 -0.02 -0.10  0.02  0.00 -0.09     0  0.11  0.04  0.13  0.09  0.08

# apply diff data on ARIMA(1,0,1)
sarima (xdiff,1,0,1)

initial  value -2.218917 
iter   2 value -2.253118
iter   3 value -2.263750
iter   4 value -2.272144
iter   5 value -2.282786
iter   6 value -2.296777
iter   7 value -2.297062
iter   8 value -2.297253
iter   9 value -2.297389
iter  10 value -2.297405
iter  11 value -2.297413
iter  12 value -2.297413
iter  13 value -2.297414
iter  13 value -2.297414
iter  13 value -2.297414
final  value -2.297414 
converged
initial  value -2.305504 
iter   2 value -2.305800
iter   3 value -2.305821
iter   4 value -2.306655
iter   5 value -2.306875
iter   6 value -2.306950
iter   7 value -2.306955
iter   8 value -2.306955
iter   8 value -2.306955
final  value -2.306955 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
         ar1      ma1   xmean
      0.3549  -0.7663  0.0072
s.e.  0.1314   0.0874  0.0032

sigma^2 estimated as 0.009885:  log likelihood = 119.88,  aic = -231.76

$degrees_of_freedom
[1] 132

$ttable
      Estimate     SE t.value p.value
ar1     0.3549 0.1314  2.7008  0.0078
ma1    -0.7663 0.0874 -8.7701  0.0000
xmean   0.0072 0.0032  2.2738  0.0246

$AIC
[1] -1.716773

$AICc
[1] -1.715416

$BIC
[1] -1.630691

sarima (xdiff,3,0,0)

initial  value -2.215090 
iter   2 value -2.289035
iter   3 value -2.306884
iter   4 value -2.308838
iter   5 value -2.309367
iter   6 value -2.309746
iter   7 value -2.309749
iter   7 value -2.309749
iter   7 value -2.309749
final  value -2.309749 
converged
initial  value -2.314672 
iter   2 value -2.314677
iter   3 value -2.314679
iter   4 value -2.314682
iter   4 value -2.314682
iter   4 value -2.314682
final  value -2.314682 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
          ar1      ar2      ar3   xmean
      -0.3669  -0.3421  -0.2363  0.0071
s.e.   0.0842   0.0849   0.0838  0.0044

sigma^2 estimated as 0.009733:  log likelihood = 120.93,  aic = -231.85

$degrees_of_freedom
[1] 131

$ttable
      Estimate     SE t.value p.value
ar1    -0.3669 0.0842 -4.3592  0.0000
ar2    -0.3421 0.0849 -4.0319  0.0001
ar3    -0.2363 0.0838 -2.8186  0.0056
xmean   0.0071 0.0044  1.6168  0.1083

$AIC
[1] -1.717413

$AICc
[1] -1.715133

$BIC
[1] -1.60981

sarima (xdiff,0,0,2)

initial  value -2.220513 
iter   2 value -2.294887
iter   3 value -2.307682
iter   4 value -2.309170
iter   5 value -2.310360
iter   6 value -2.311251
iter   7 value -2.311636
iter   8 value -2.311648
iter   9 value -2.311649
iter   9 value -2.311649
iter   9 value -2.311649
final  value -2.311649 
converged
initial  value -2.310187 
iter   2 value -2.310197
iter   3 value -2.310199
iter   4 value -2.310201
iter   5 value -2.310202
iter   5 value -2.310202
iter   5 value -2.310202
final  value -2.310202 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
          ma1      ma2   xmean
      -0.3984  -0.2173  0.0072
s.e.   0.0808   0.0768  0.0033

sigma^2 estimated as 0.00982:  log likelihood = 120.32,  aic = -232.64

$degrees_of_freedom
[1] 132

$ttable
      Estimate     SE t.value p.value
ma1    -0.3984 0.0808 -4.9313  0.0000
ma2    -0.2173 0.0768 -2.8303  0.0054
xmean   0.0072 0.0033  2.1463  0.0337

$AIC
[1] -1.723268

$AICc
[1] -1.721911

$BIC
[1] -1.637185

Symperasmata Apo tin efarmogi ton 3wn parapano montelon prokiptei oti ola kanoun kalo fit dioti: - kai ta 3 exoun ACF residuals Mi statistika simantika - kai Ljung-Box me statistika simantika p-values.

Wstoso gia provlepsi epilegoume to montelo ARIMA(1,1,1) prokimenou na pragmatopoiisoume provlepsi gia ta epomena 10 xronia dioti exei tous xamiloterous AIC kai BIC diktes

sarima.for(globtemp, 10, 1, 1, 1)

$pred
Time Series:
Start = 2016 
End = 2025 
Frequency = 1 
 [1] 0.8031838 0.7841177 0.7819961 0.7858872 0.7919120 0.7986940 0.8057447 0.8128907 0.8200705 0.8272623

$se
Time Series:
Start = 2016 
End = 2025 
Frequency = 1 
 [1] 0.09942219 0.11536298 0.12346802 0.12943468 0.13462409 0.13945017 0.14405681 0.14850078 0.15280874 0.15699617

Erotima 5)

Θεωρείστε τα δεδομένα chicken Προσαρμόστε αυτά σε κατάλληλο μοντέλο ARIMA και κάντε πρόβλεψη για τους επόμενους 12 μήνες.

Analisi Dedomenon Poultry (chicken), Whole bird spot price, Georgia docks, US cents per pound

# load globtemp dataset and make sure that it is a time series 
x=ts(chicken)

# plot oil data
plot.ts(chicken, main="Chicken Dataset")

# plot oil lag 1
lag1.plot(chicken,1)

# plot ACF and PACF 
acf2(x)

     [,1]  [,2]  [,3]  [,4]  [,5]  [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  0.99  0.97  0.95  0.93  0.91  0.89 0.87 0.86 0.84  0.82  0.80  0.78  0.75  0.73  0.71  0.68  0.66  0.63  0.61  0.59  0.57  0.55  0.53   0.5
PACF 0.99 -0.11 -0.08 -0.03 -0.02 -0.01 0.00 0.00 0.00 -0.01 -0.05 -0.06 -0.07 -0.04 -0.03  0.00  0.00  0.02  0.01  0.02 -0.02  0.00 -0.01   0.0

Differecing

Paratiroume oti ta dedomena exoun auksitiki tasi ta dedomena opote tha efarmosoume 1es diafores. Apo ta ACF kai PACF paratiroume oti to ACF akolouthi fthinei ekthetika akolouthontas imitonoidi motivo kai episis paratiroume 2 spikes sto lag 1 kai 2 gia to PACF. Endiksis pou mas odigoun na xrisimopoiisoume ena AR(2) montelo pano sta dedomena pou proekipsan apo tis diafores

Symperasmata Parolo pou ta ACF kai PACF mas odigoun stin efarmogi enos AR(2) montelou gia tin proti diafora. Ta apotelesmata den einai poli entharintika dioti: - Paratiroune orismena statistika simanika residuals - Kathos kai mi statistika P-values to ljung-Box plot.

xdiff=diff(x,1)
# plot globtemp data
plot.ts(xdiff, main="Chicken Dataset After First Diffs")

# plot ACF and PACF for first diff data
acf2(xdiff)

     [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24]
ACF  0.72  0.39  0.09 -0.07 -0.16 -0.20 -0.27 -0.23 -0.11  0.09  0.26  0.33  0.20  0.07 -0.03 -0.10 -0.19 -0.25 -0.29 -0.20 -0.08  0.08  0.16  0.18
PACF 0.72 -0.29 -0.14  0.03 -0.10 -0.06 -0.19  0.12  0.10  0.16  0.09  0.00 -0.22  0.03  0.03 -0.11 -0.09  0.01 -0.03  0.07 -0.04  0.06 -0.05  0.02

# apply diff data on ARIMA(2,0,0)
sarima (xdiff,2,0,0)

initial  value 0.001863 
iter   2 value -0.156034
iter   3 value -0.359181
iter   4 value -0.424164
iter   5 value -0.430212
iter   6 value -0.432744
iter   7 value -0.432747
iter   8 value -0.432749
iter   9 value -0.432749
iter  10 value -0.432751
iter  11 value -0.432752
iter  12 value -0.432752
iter  13 value -0.432752
iter  13 value -0.432752
iter  13 value -0.432752
final  value -0.432752 
converged
initial  value -0.420882 
iter   2 value -0.420934
iter   3 value -0.420935
iter   4 value -0.420936
iter   5 value -0.420936
iter   6 value -0.420937
iter   6 value -0.420937
iter   6 value -0.420937
final  value -0.420937 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
         ar1      ar2   xmean
      0.9494  -0.3069  0.2632
s.e.  0.0717   0.0718  0.1362

sigma^2 estimated as 0.4286:  log likelihood = -178.64,  aic = 365.28

$degrees_of_freedom
[1] 176

$ttable
      Estimate     SE t.value p.value
ar1     0.9494 0.0717 13.2339  0.0000
ar2    -0.3069 0.0718 -4.2723  0.0000
xmean   0.2632 0.1362  1.9328  0.0549

$AIC
[1] 2.040697

$AICc
[1] 2.041463

$BIC
[1] 2.111923

Aksiologisi montelou To montelo AR(2) Den douleuei opos tha itan anamenomeno kai auto giati paratirite oti sta dedomena mas meta tin efarmogi ton proton diaforon iparxei epoxikotita kai poio sigkekrimena monthlly seasonality. Auto fenete kai apo to grafimia tis xronosiras twn prwton diaforwn alla kai apo to ACF grafima to opoio akolouhti imitonoidi morfi me periodo ena lag.

Gia tin eksalipsi autis tis epoxikotitas tha xrisimopoiisoume ena seasonal arima SARIMA(2,0,0)x(1,0,0) me S = 12.

To opoio fenete na douleuei poli kala me ta dedomena mas kathos: - Ola ta residuals einai mi statistika simantika - kai ola ta p-values sto Ljung-Box plot einai statistika simantika.

# apply diff data on ARIMA(2,0,0)
sarima (xdiff,2,0,0, P=1,D=0,Q=0, S=12)

initial  value 0.015039 
iter   2 value -0.226398
iter   3 value -0.412955
iter   4 value -0.460882
iter   5 value -0.470787
iter   6 value -0.471082
iter   7 value -0.471088
iter   8 value -0.471090
iter   9 value -0.471092
iter  10 value -0.471095
iter  11 value -0.471095
iter  12 value -0.471096
iter  13 value -0.471096
iter  14 value -0.471096
iter  15 value -0.471097
iter  16 value -0.471097
iter  16 value -0.471097
iter  16 value -0.471097
final  value -0.471097 
converged
initial  value -0.473584 
iter   2 value -0.473663
iter   3 value -0.473721
iter   4 value -0.473822
iter   5 value -0.473871
iter   6 value -0.473884
iter   7 value -0.473885
iter   8 value -0.473885
iter   8 value -0.473885
iter   8 value -0.473885
final  value -0.473885 
converged
$fit

Call:
arima(x = xdata, order = c(p, d, q), seasonal = list(order = c(P, D, Q), period = S), 
    xreg = xmean, include.mean = FALSE, transform.pars = trans, fixed = fixed, 
    optim.control = list(trace = trc, REPORT = 1, reltol = tol))

Coefficients:
         ar1      ar2    sar1   xmean
      0.9154  -0.2494  0.3237  0.2353
s.e.  0.0733   0.0739  0.0715  0.1973

sigma^2 estimated as 0.3828:  log likelihood = -169.16,  aic = 348.33

$degrees_of_freedom
[1] 175

$ttable
      Estimate     SE t.value p.value
ar1     0.9154 0.0733 12.4955  0.0000
ar2    -0.2494 0.0739 -3.3728  0.0009
sar1    0.3237 0.0715  4.5238  0.0000
xmean   0.2353 0.1973  1.1924  0.2347

$AIC
[1] 1.945973

$AICc
[1] 1.947257

$BIC
[1] 2.035006